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ABSTRACT 


A  primitive  equation  model  is  used  to  simulate  mesoscale  ocean  variability 
to  verify  the  digital  filter  initialization  (DPI)  method  of  Lynch  and  Huang.  The 
model  is  initialized  with  two  different  analytic  density  fields,  a  linear  Rossby  wave 
and  a  baroclinically  unstable  flow,  and  then  integrated  forward  to  produce  control 
runs.  Both  simulations  resulted  in  measurable  ageostrophic  currents  and  vertical 
motions.  The  density  field  at  the  end  of  the  control  runs  was  extracted  and  used 
by  DPI  to  regenerate  the  control  currents.  Combinations  of  different  DPI 
integration  times  and  repeated  DPI  passes,  were  used.  The  normalized  rms  error 
between  the  vertical  velocities  from  the  control  runs  and  from  DPI,  showed  that 
the  12-hr  and  3  pass  combination  had  the  greatest  accuracy.  The  normalized  error 
was  less  than  0.15,  except  near  the  bottom  of  the  domain  and  also  at  the  surface 
for  the  shallow  baroclinically  unstable  flow.  This  was  attributed  to  the  neglect  of 
friction  in  the  DPI  processes,  and  to  somewhat  poor  vertical  resolution  near  the 
surface  for  the  unstable  case.  The  errors  were  small  enough  to  be  confident  that 
DPI  accurately  recovered  the  adiabatic,  frictionless  part  of  the  control’s  currents. 
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I .  IMTBODUCTZON 


In  recent  years  the  U.S.  Navy  has  shifted  its  attention 
and  operation  from  open-ocean  warfare,  to  expeditionary  roles 
in  the  "littoral"  environment  of  the  coastlines  of  the  world. 
The  new  emphasis  was  first  set  forth  in  the  U.S.  Navy's  white 
paper,  "...From  the  Sea,"  (O'Keefe  et  al.  1992)  and  amplified 
in  the  follow-on  white  paper,  "Forward. .  .From  the  Sea"  (Dalton 
et  al.  1994).  The  needs  of  a  naval  expeditionary  force  are 
different  from  those  of  a  deep  water  force.  Anti-Submarine 
warfare  becomes  a  difficult  challenge  in  the  coastal  waters, 
where  sound  propagation  is  very  complex  and  not  modeled  very 
well.  Mine  warfare  becomes  increasingly  important  with  a 
desire  for  an  increased  ability  to  detect  moored  mines,  and  a 
complex  environment  that  poses  unique  challenges  for  all 
forces  that  operate  within  these  areas.  This  focus  on  the 
littoral  has  caused  increased  attention  on  numerical  coastal 
prediction  models  with  a  desire  to  increase  the  accuracy  of 
their  forecasts. 

One  direct  method  to  raise  the  accuracy  of  the  product  is 
to  increase  the  accuracy  of  the  initial  conditions.  This 
goal  could  be  satisfied  by  the  digital  filter  initialization 
(DFI)  method  of  Lynch  and  Huang  (1992).  DFI  would  accomplish 
this  by  reducing  the  gravity  wave  (high  frequency)  noise  at 
the  start  of  a  model  forecast  cycle.  Gravity  wave  noise  is  a 
result  of  the  imbalances  in  the  horizontal  equations  of 
motion.  By  removing  these  imbalances  frcMn  the  model 
additional  shoc)ca,  from  the  imbalances  between  the  model 
dynamics  and  the  vertical  circulation,  could  be  eliminated. 
The  initialization  procedure  must  also  emulate  the  effects  of 
the  geostrophic  adjustment  process,  so  that  the  initialized 
currents  have  a  divergent  conponent ,  resulting  in  vertical 
motions. 

Horizontal  velocities  can  be  measured  directly  with 
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current  meters,  drifting  buoys,  etc.  However,  the  accuracy  of 
these  measurements  is  not  sufficient  to  directly  determine  the 
vertical  velocity  through  mass  continuity.  The  difficulty 
with  accurately  initializing  a  3-dimensional  circulation  model 
therefore  lies  with  how  the  vertical  velocities  are 
determined.  These  must  be  indirectly  inferred  or  diagnosed 
with  the  accuracy  depending  upon  the  assuirptions,  methods 
used,  and  the  manner  the  data  were  collected. 

Many  operational  weather  prediction  centers,  including 
the  Fleet  Numerical  Oceanography  Center  (FNMOC),  use  nonlinear 
normal  mode  initialization  (NNMI)  to  begin  their  forecast 
models.  However,  Lynch  and  Huang  (1992)  recently  showed  that 
DFI  is  equivalent  (if  not  superior)  to  NNMI  and  is  far  easier 
to  apply. 

The  purpose  of  this  thesis  is  to  validate  the  DFI  method 
within  the  context  of  a  coastal  ocean  prediction  model.  If 
the  DFI  method  is  quantitatively  validated,  it  can  be  used  as 
an  initialization  method  for  numerical  ocean  prediction 
models.  It  can  also  be  utilized  to  diagnose  the  vertical 
circulation  and  ageostrophic  velocity  from  observed  data  sets 
and  thereby  used  in  dynamical  studies. 

A.  BACKSROCniD 

1.  Nathods  to  Detoxaine  Dlverganoe 

Several  methods  have  been  employed  to  diagnose  the 
divergent  component  of  the  current  and  its  associated  vertical 
velocity.  The  following  is  a  brief  review  of  a  selection  of 
several  methods  used  to  reveal  the  vertical  velocity 
structure. 

a.  Umthodm  Bmmmd  on  Dynmmxcml  and  TbmxwodynamxaMl 

EquMtlooM 

In  his  paper  Stress  (1994)  reviewed  a  set  of 
nvimerical  methods.  He  analyzed  and  ranked  the  isopycnic 
advection  equation,  the  vorticity  advection  equation  and  the 
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Q-vector  version  of  the  omega  equation,  as  they  were  applied 
to  estimate  the  mesoscale  distribution  of  vertical  motions. 
Strass  used  the  mesoscale  term  to  refer  to  scales  as  small  as 
the  Rossby  radius  of  deformation,  which  is  typically  on  the 
order  of  25  km.  The  methods  were  tested  to  determine  how 
their  computed  vertical  velocities  compared  to  "observed",  and 
if  the  w  velocity  was  consistent  with  the  model  variables. 
"Observed"  values  were  determined  from  controlled  numerical 
solutions.  The  results  showed  that  the  Q-vector  method  was 
best  for  diagnosis  of  the  horizontal  pattern  of  vertical 
motion.  Nevertheless,  when  the  integration  time  advanced 
causing  the  vertical  velocity  to  have  its  strongest  variance, 
the  coherence  decreased  for  wavelengths  less  than  10  km  {high 
wavenumbers)  (Strass  1994) . 

Advection  of  vorticity  was  found  to  be  the  second 
best  method  for  predicting  the  model's  vertical  velocity. 
However,  here  too,  was  a  drop  of  coherence,  at  even  higher 
wavenumbers  than  those  found  in  the  Q— vector  method, 
(beginning  at  30  km) .  Additionally  a  looser  relationship  to 
the  modeled  w,  was  found  concerning  the  phase. 

Strass  (1994)  found  the  least  accurate  model  was  the 
isopycnic  advection  model,  which  diagnoses  w  from  a  steady 
state  thermodynamic  equation  (assuming  the  horizontal  and 
vertical  advection  balanced) .  It  reproduced  the  large  scale 
distribution  of  vertical  motion,  but  it  was  deficient  in  the 
smaller  scales,  as  no  individual  upwelling  or  downwelling 
areas  were  correctly  identified.  In  fact  there  was  a 
progressive  phase  error  over  time  that  settled  near  180 
degrees  out  of  phase  (i.e.,  upwelling  predicted  for  actual 
downwelling  and  vice  versa) .  The  overall  correlation  was 
therefore  lower  than  the  other  two  methods.  The  stationarity 
assumption  becomes  fatal  for  the  oceanic  synoptic  scale. 

Vertical  velocity  distribution  patterns  were  best 
predicted  by  the  Q-vectors  of  the  quasi-geostrophic  (QG)  omega 
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equation,  and  the  amplitude  had  an  error  of  only  a  few  tens  of 
percent.  However,  to  convert  the  divergence  from  the 
voriticity  advection,  or  the  convergence  from  the  Q-vector,  or 
the  Q-vector  potential  into  vertical  motion,  a  scaling  factor 
had  to  be  used  (Stress  1994).  This  result  was  expected  since 
the  Q-vector  form  of  the  QG  omega  equation  is  known  to  be  the 
more  accurate  method  (Holton  1992) .  The  QG  canega  equation  is 
derived  from  consistent  vorticity  and  thermodynamic  equations, 
with  no  assumptions  of  steady  state.  It  is  however,  only 
valid  for  QG  motions,  i.e.,  those  scales  of  motions  for  which 
the  Rossby  number  is  small. 

b.  Drifting  Baoym 

Several  studies  have  used  drifting  buoys  to  directly 
find  divergence  and  indirectly  imply  vertical  motion.  Swenson 
et  al.  (1988)  used  mixed  layer  drifters'  trajectories  and 
temperature  measurements  to  arrive  at  the  vertical  velocity. 
They  did  not  estimate  divergence  directly  from  the  horizontal 
con^nent  of  the  flow,  but  rather  indirectly  from  the 
vorticity  budget.  This  required  accurate  estimations  of  the 
relative  vorticity,  which  allowed  the  time  derivatives  of  the 
buoys  motion  to  be  calculated  (Swenson  et  al.  1988).  The 
short  comings  of  this  were:  that  the  estimated  relative 
vorticity  was  taken  from  a  subjectively  chosen  cluster  of 
buoys;  and  that  analyzing  their  relative  motions  and  resulting 
vorticity  changes  were  done  for  only  small  scales  (12-15  km) . 

Lagrangian  drifters,  satellite  images,  direct 
sampling  and  hydrographic  surveys  were  used  by  Paduan  and 
Niiler  (1990),  to  track  water  parcels  at  a  given  depth  in 
California  coastal  waters.  From  their  movement,  estimates  of 
the  vertical  velocity  were  made.  "To  attach  physical 
relevance  to  the  velocity  gradient  estimated  using  the  cluster 
method,  required  making  the  assumption  that  there  exists  a 
separation  in  scales  between  cluster-scale  and  sub-cluster 
scale  motions"  (Paduan,  Niiler  1990) .  They  used  the  shallow 
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water  approximation,  and  assumed  that  the  wind  driven  curl 
effects  were  not  important  in  the  well-mixed  layer. 
Nevertheless,  they  did  allow  Ekman  transport.  The  direct 
estimates  of  divergence  were  unreliable  and  suggested  a 
measure  of  uncertainty  in  the  cluster  method,  when  applied  to 
determining  divergence.  This  uncertainty  exists  particularly 
when  small  numbers  of  drifters  are  used,  and  when  they  spanned 
a  maximum  in  the  flow  field.  Relative  vorticity  was  better 
estimated  and  used  (with  the  vorticity  budget)  to  locate  areas 
of  upwelling  and  downwelling. 

Isopycnal  RAFOS  floats  were  used  off  Cape  Hatteras 
by  Bower  (1989)  to  examine  convergence  and  divergence.  She 
estimated  the  horizontal  divergence  from  the  average  rate  of 
change  of  the  absolute  vorticity  with  time  between  the 
extremes  of  Gulf  Stream  meanders  (Bower  1989) .  The  change  in 
time  (At)  was  the  time  it  took  the  float  to  traverse  the 
segment  over  which  the  absolute  vorticity  was  calculated. 
Uncertainties  in  the  estimations  of  the  divergence  were  a 
function  of  At  and  absolute  vorticity.  Float  trajectories 
were  used  to  estimate  the  sign  of  the  horizontal  divergence 
and  the  magnitudes  of  that  divergence  from  the  vorticity 
budget.  The  results  were  questionable  due  to  uncertainties  in 
estimating  the  vorticity  associated  with  the  lateral  shear 
along  a  density  surface.  Thus  they  were  only  able  to  give 
upper  and  lower  limits  of  divergence  in  the  main  thermocline 
(Bower  1989) . 

o.  Txmakinff  thitxxmatm  mad.  CStloxe^pbyll 

Beside  the  use  of  drifting  buoys,  the  employment  of 
satellite  imagery,  and  san^ling  of  the  water  have  been  applied 
to  track  nutrients,  chlorophyll,  etc.,  to  estimate  the 
vertical  transportation  in  the  ocean  (Kadko  et  al.  1991).  By 
tracking  an  isotope  related  to  chlorophyll,  Kadko  et  al. 
(1991)  tried  to  determine  how  deep  the  chlorophyll  would  sink. 
They  looked  at  the  depth  chlorophyll  would  sink  to  the  half— 
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life  of  the  isotope  found  at  that  depth.  The  isotope  was 
assumed  to  combine  with  the  chlorophyll  at  the  surface  and 
with  its  rapid  half-life  would  only  exist  a  short  time  as  the 
chlorophyll  began  to  sink.  Chlorophyll  was  found  at  depths 
deeper  than  could  be  reached  by  the  calculated  sink  rate.  The 
presence  of  the  isotope  also  showed  that  the  chlorophyll  was 
not  being  created  by  the  phytoplankton  at  that  depth.  Kadko 
et  al.  (1991)  concluded  it  must  have  been  the  result  of  rapid 
vertical  transport  with  values  estimated  to  be  5-20  m  d“^  in 
the  California  current  eddies. 

Washburn  et  al.  (1991)  traced  high  concentrates  of 
phytoplankton,  and  found  phytoplankton  at  depths  without 
enough  light  to  grow  and  reproduce.  They  used  the 
distribution  of  chlorophyll  fluorescence,  to  give  vertical 
resolution.  Washburn  et  al.  (1991)  concluded  that  particle 
sinking  rates  were  not  the  dominate  process  nor  was  the 
resuspension  of  bottom  sediments  that  included  phytoplankton. 
They  resolved  that  to  get  chlorophyll  at  the  observed  depth, 
there  had  to  be  a  subduction  of  surface  waters.  Washburn  et 
al.  (1991)  also  suggested,  however  that  some  vertical  movement 
could  be  the  result  of  mixing  along  sloping  isopycnals.  The 
vertical  displacement  based  on  the  geostrophic  flow  on  sloping 
density  surfaces  underestimated  the  vertical  motions.  This 
underestimation  resulted  from  the  calculations  not  going  to 
the  nearshore,  where  the  isopycnals  had  the  largest  slope  into 
the  higher  velocity  nearshore  flow  (Washburn  et  al.  1991). 

d.  DFI  in  tbm  Atmo^pbmxm 

Lynch  and  Huang  (1992),  used  the  digital  filter 
initialization  (DFI)  method  to  initialize  the  High  Resolution 
Limited  Area  Model  (HRLAM) .  The  initialized  fields  represent 
dynamically  balanced  model  variables  with  a  realistic 
distribution  of  the  ageostrophic  flow  and  vertical  velocity. 
HRLAM  is  a  primitive  equation  forecast  model  of  the 
atmosphere.  A  lowpass  digital  filter  was  applied  to  the  time 
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series  generated  by  a  simple  3-hour  forward  and  3-hour 
backward  model  integration  (6-hour  total  time  span) .  This 
lowpass  filtered  time  series,  centered  at  the  initial  time, 
represented  the  dynamically  balanced  (initialized)  fields  of 
the  model.  Lynch  and  Huang  (1992)  also  showed  that  this 
balanced  field  is  equivalent  to  that  obtained  by  the  more 
complicated  NNMI.  A  doubling  of  the  total  time  span  to  12 
hours  led  to  a  dandling  of  the  meteorological  modes.  The 
forward  and  backward  integrations  were  accomplished  with  the 
physics  turned  off,  to  allow  no  irreversible  processes  to  be 
integrated  backwards.  The  boundaries  were  held  constant,  but 
the  fields  inside  were  allowed  to  vary.  Lynch  and  Huang 
(1992)  showed  that  DFI  eliminated  initial  gravity  wave  noise 
from  the  models  used  in  atmospheric  forecasting. 

B.  (»JBCTIVB 

The  objective  of  this  study  is  to  prove  that  DFI  can 
diagnose  accurately  the  secondary  circulations  associated  with 
quasi-geostrophic,  and  other  higher  order,  dynamical  balances. 
By  illustrating  the  accuracy  of  the  DFI  method,  it  can  then  be 
applied  with  confidence,  to  better  understand  and  predict 
mesoscale  features  in  the  littoral  regions  of  the  oceans. 


C.  MOOBL  MID  MBTBOOS 

Following  the  approach  of  Strass  (1994),  the  method  used 
to  test  DFI  was:  to  extract  the  density  field  from  a 
controlled  numerical  model  simulation;  apply  the  DFI  method  to 
initialize  the  currents  (including  the  vertical  velocity);  and 
finally  verify  the  resulting  initialized  currents  against  the 
known  currents  from  the  control  simulation.  Two  different 
simulations  are  carried  out.  One  corresponded  to  linear 
Rossby  waves  and  the  other  corresponded  to  a  field  of  growing 
baroclinic  disturbances  typical  of  the  coastal  ocean. 
Starting  a  primitive  equation  numerical  ocean  model  from  an 
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initial  condition  of  rest  or  geostrophic  balance,  causes  the 
currents  and  mass  fields  to  undergo  relatively  high  frequency 
oscillations,  as  a  part  of  the  geostrophic  adjustment  process. 
When  the  Initial  condition  is  at  a  state  of  rest,  the  model 
has  many  more  oscillations  than  when  it  is  started  from  the 
geostrophic  balance,  and  the  adjustment  process  causes  the 
geostrophic  currents  to  be  spun  up.  The  ageostrophic  currents 
are  created  by  the  geostrophic  adjustment  process  when 
starting  from  the  geostrophic  balance  case.  The 
geostrophic/ageostrophic  currents  and  the  secondary 
circulation  make  up  the  "slow  manifold"  (Lorenz  1992)  in  the 
model.  The  "slow  manifold"  is  that  part  of  the  oceanic  state 
that  may  be  considered  determinably  predictable.  Lynch  and 
Huang  (1992)  showed  that  a  lowpass  digital  filter  applied  to 
the  model's  variables  removed  the  high  frequency  oscillation 
(i.e.,  gravity  waves),  and  returned  the  "slow  manifold"  at  the 
analysis  time.  These  variables  were  generated  by  short  term 
forward  and  backward  integrations  that  start  from  an 
uninitialized  analysis. 

The  numerical  model  used  is  a  multi-level,  hydrostatic, 
nonlinear  primitive  equation  model.  The  vertically  integrated 
currents  are  nondivergent  (i.e.,  "rigid  lid").  There  are  20 
levels  in  the  vertical,  however  the  levels  where  w  is 
calculated  are  staggered  between  the  levels  where  u  and  v  are 
calculated.  For  a  given  level  the  location  of  w  is  half  a 
level  lower  than  that  for  u  and  v.  In  the  horizontal  a  4  km 
grid  spacing  was  used.  The  model  is  a  cyclic  zonal  channel 
with  insulated  free-slip  walls  at  the  north  and  south 
boundaries,  and  a  flat  bottom.  The  model  domain  is  224  km 
from  the  east  to  the  west  and  224  km  frcan  the  south  to  the 
north,  with  a  total  depth  of  4000  m. 

The  DFI  method  used  is  similar  to  Lynch  and  Huang's 
(1992)  approach.  A  time  series  of  the  model  variables  is 
generated  by  a  short  term  forward  and  backward  integration 
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from  the  initial  time  (Figure  1).  The  integrations  are  done 
without  any  forcing  or  friction.  The  resulting  time-series  is 
of  length  2T,  and  a  lowpass  filter  is  then  applied  to  it. 


t  =  -  r  t  =  0  t  =  T 


Figure  1.  Schematic  representation  of  the  DFI  method  developed 
by  Lynch  and  Huang  1992.  Starting  from  a  given  initial 
condition,  the  model  is  integrated  forward  and  backward. 

Thus,  if  q  represents  any  one  of  the  model's  variables, 
then  the  lowpass  filtered  value,  corresponding  to  the 
dynamically  balanced  state  at  the  time  t  =  0,  is  denoted  by 
qj,  which  given  by: 

g,  =  E 

n~-N 

Where  h„  is  the  lowpass  filter  coefficient  (Lynch  and  Huang 
1992),  and  q„  =  q(nAt)  is  the  discrete  value  of  the  model 
variable  at  time  t  =  nAt,  for  — T  s  t  s  T  and  T  —  NAt. 

Lynch  and  Huang  (1992)  showed  that  the  integration  and 
lowpass  filter  removed  the  high  frequency  noise  and  produced 
dynamically  balanced  model  fields.  These  fields  now  have  very 
little  noise  associated  with  them,  when  used  to  initialize  the 
model  and  in  the  forward  integration  of  the  model's  forecast. 
In  this  study  the  accuracy  with  which  the  DFI  can  diagnose  the 
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known  currents,  by  starting  only  with  the  known  density  field, 
is  tested.  In  addition,  the  optimum  value  of  the  filter  span 
T  is  determined  and  the  usefulness  of  repeated  initialization 
is  also  examined.  A  repeated  initialization  is  one  where  the 
initialized  variables,  q„  are  used  as  initial  conditions  for 
subsequent  forward  and  backward  integration  plus  filter 
application.  In  fact  several  such  repeated  initialization 
steps  were  used  to  examine  the  resulting  convergence  to  the 
known  solution  (control  currents) . 
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IX .  RESULTS 


A.  CONTROL  RUN  -  LINEAR  ROSSBT  NAVES 

1 .  Seti9 

The  model  was  first  initialized  using  DFI  with  a  12-hour 
forward  and  a  12-hour  backward  integration  {2T  -  24  hr  total 
time  span)  and  it  was  repeated  3  times.  The  DFI  began  at  t  “ 
0  with  a  large  scale,  small  amplitude  density  disturbance  with 
currents  that  were  in  geostrophic  and  hydrostatic  balance. 
This  density  disturbance  was  defined  in  the  model  by  the 
following  equation: 


xm 


dist  ■  (sin(6*E* 


))♦ 


xm 

(1  +  COS  (2*e* - )  $ 

ym 


(2) 


where  xm  and  ym  are  both  equal  to  224  km,  and  are  the  x-  and 
y-  scales  of  the  domain.  It  has  a  sine  wave  pattern  with  a 
wavenumber  of  3  in  the  x-direction.  Along  the  y-direction  the 
density  disturbance  has  a  cosine  pattern  of  wavenumber  1. 
This  results  in  dist  -  0  at  the  north  (y  -  ym)  and  south  (y  - 
0)  boundaries,  and  a  maximum  in  the  center  of  the  y-axis  (y  = 
ym/2)  .  The  disturbance  in  Equation  (2)  was  then  multiplied  by 
a  "form  function"  to  give  it  a  vertical  profile  that  was 
constant  from  the  surface  to  a  depth  of  37.5  m.  Following 
that  it  underwent  a  smooth  transition  to  an  exponential 
profile  beginning  at  75  m.  The  exponential  profile  had  a 
scale  depth  of  240  m. 

This  density  disturbance  was  superimposed  on  a  stratified 
zonal  shear  flow  that  was  also  in  geostrophic  and  hydrostatic 
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balance.  That  current's  density  field  was  defined  by  the 
hyperbolic  tangent  of  (y/yl) ,  where  yl  *  ym*10^.  A  nearly 
constant  positive  northward  density  gradient  was  created,  such 
that  the  northern  boundary  of  the  model  domain  had  a  higher 
density  than  the  southern  boundary.  The  horizontal  structure 
was  also  multiplied  by  the  same  z-dependent  function  as  the 
density  disturbance  (constant  in  z  above  37.5  m,  and  then  a 
transition  to  an  exponential  function  at  75  m,  with  a  240  m 
scale  depth) .  The  density  gradient  produced  a  nearly  uniform 
zonal  shear  flow  (with  a  corresponding  exponential  profile  in 
z)  that  was  centered  in  the  middle  of  the  y-axis. 

The  density  disturbance  and  zonal  flow  were  used  in  the 
DFI  with  a  12-hr  integration  time  span  and  3  repeated  passes. 
From  the  DFI  the  density  field  and  the  associated  currents 
were  used  to  initialized  the  model.  After  the  initialization, 
the  model  was  integrated  forward  for  10  days,  which  allowed 
the  large  scale  disturbance  to  evolve  and  propagate.  The 
disturbance  produced  by  Equation  (2),  evolved  as  a  quasi- 
linear  first  baroclinic  mode  Rossby  wave  that  was 
representative  of  a  typical  oceanic  "slow  mode".  The  model  is 
on  a  f-plane  and  therefore,  the  "beta"  effect  was  due  to  the 
wea)c  northward  gradient  of  the  mean  relative  vorticity.  The 
instantaneous  density  field  was  then  extracted  from  the  10th 
day  of  the  10-day  control  integration.  The  day  10  density 
field  was  then  used  with  the  DFI  procedure  with  the  purpose  of 
trying  to  recover  the  original  day  10  (slightly 
nongeostrophic)  control  currents.  One  set  of  DFI  tests 
specified  the  initial  horizontal  velocities  to  be  the 
geos trophic  velocity,  calculated  from  the  day  10  control 
density  structure.  For  a  second  set  of  tests,  the  initial 
horizontal  velocity  was  zero.  The  DFI  was  performed  with  four 
different  integration  times  (T  =  6,  9,  12,  and  15  hours),  and 
up  to  three  repeated  initialization  passes  for  each 
integration  time. 
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2 .  Control  Slnulation 

The  linear  Rossby  wave  simulation  resulted  in  the 
distribution  of  divergence,  vorticity  and  ageostrophic 
currents  shown  in  Figure  2.  In  the  upper  ocean,  the  maximum 
upward  vertical  velocity  is  co-located  on  the  maximum 
divergence,  and  the  maximum  downward  velocity  is  co-located  on 
the  maximum  convergence  (Figure  2. a,  and  3. a).  The  vector 
form  of  the  ageostrophic  currents  illustrates  this  convergence 
and  divergence  (Figure  2.c) .  A  regular  distribution  was  also 
found  for  the  u  and  v  velocities,  with  a  slight  deformation 
near  the  north  and  south  boundaries  (Figures  3.b,  and  3.c). 
The  tilted  vertical  structure  of  w,  produced  by  the 
disturbance,  is  seen  in  Figure  4 .  The  tilt  of  the  w  and  v 
(not  shown)  vertical  structures  in  the  direction  of  the  mean 
flow  suggest  a  possible  dan^ening  of  the  fields  (Pedlosky 
1964) . 

The  secondary  circulation  associated  with  the  Rossby  wave 
disturbance  is  in  quasi-geostrophic  balance.  There  is 
horizontal  convergence,  and  horizontal  divergence.  The 
horizontal  convergence  has  an  associated  pattern  of 
downwelling  that  is  downstream  of  the  ridges,  while  the 
horizontal  divergence  has  its  associated  upwelling  pattern 
located  downstream  of  the  troughs.  The  focus  of  this  paper  is 
on  the  divergent  part  of  the  secondary  circulation  and  its 
associated  vertical  velocity.  These  represent  a  large  part  of 
the  ageostrophic  current  signal,  and  are  an  important  part  of 
ocean  dynamics. 

3 .  DFX  V«rlf icatlon 

As  described  above,  the  DFI  method  was  used  to  try  to 
diagnose  the  currents  at  day  10  from  the  control  run.  This  is 
done  by  only  using  the  (control)  density  field  at  that  time. 
The  accuracy  of  the  DFI  solution,  as  a  function  of  the 
different  integration  time  spans,  was  first  examined  by  way  of 
a  normalized  root  mean  square  (rms)  error  (NRE) .  This 
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X  axis 

Figure  4.  Rossby  Wave.  W  vertical  profile  (m  d“^)  . 

consisted  of  the  rms  of  the  difference  between  the  vertical 
velocity  computed  from  the  initialized  (DFI)  currents,  and  the 
true  (control  at  day  10)  vertical  velocity,  which  is  then 
normalized  by  the  rms  of  the  control's  vertical  velocity.  A 
percentage  of  the  error  with  regard  to  the  true  vertical 
velocity  signal  is  returned  from  Equation  (3) . 


[  ^  . 


In  this  equation  the  w,  represents  the  filtered  (DFI)  vertical 
velocity,  w  is  the  control's  vertical  velocity,  and  the 
summation  indicates  a  horizontal  average.  In  calculating  the 
average,  the  outer  5  grid  points  in  the  x  and  y  directions 
were  not  used.  This  eliminates  from  the  verification  any 
problems  that  might  occur  on  the  boundaries  of  the  model  due 
to  the  existence  of  any  'non-slow  mode'  behavior  (gravity 
waves)  located  there. 

A  verification  was  first  made  for  each  of  the  four 
integration  times  (T  =  6,  9,  12,  and  15  hrs),  with  just  one 
DFI  pass.  For  all  integration  times,  the  control's  u  and  v 
velocity  fields  were  almost  completely  recovered  (small  NRE) 


by  the  first  DFI  pass.  However,  the  w  field  required  the 
integration  time  to  be  lengthened,  from  T  =  6  to  T  =  15  hrs, 
in  order  for  more  of  it  to  be  returned  (Figure  5) . 


Figure  5.  Normalized  RMS  Error.  First  DFI  pass  with  different 
integration  times.  Solid  line,  6-hr;  dot  dash  line,  9-hr;  dot 
line,  12-hr;  dash  line,  15-hr. 


As  shown  in  Figure  5,  there  is  a  dramatic  decrease  in  the 
normalized  rms  error  when  T  is  increased  from  6-hrs  to  9-hrs. 
From  the  9-hr  to  15-hr  integrations  there  are  progressively 
smaller  steps  in  advancing  the  accuracy.  The  non-convergence 
(NRE  ~  1)  in  the  deeper  levels  is  addressed  below. 

The  results  from  repeating  the  DFI  procedure  (repeated 
"initialization")  up  to  three  times,  were  then  compared  for 
each  integration  time.  For  the  four  different  integration 
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times  (T  -=  6,  9,  12,  and  15  hrs),  the  results  show  that  the 
normalized  rms  error  between  the  w's,  decreased  with  the 
number  of  DFI  passes  (Figures  6,  and  7).  The  least 
improvement  concerning  the  normalized  rms  error  was  found  in 
the  second  to  third  pass.  Notable  are  the  second  and  third 
passes  of  the  12-hr  and  15-hr  integration  times  (Figure  7). 
This  suggests  that  after  the  second  DFI  pass  almost  all  of  the 
control  vertical  velocity  has  been  recovered.  When  comparing 
pass  nxmber  3  for  the  different  time  spans,  the  longer 
integration  time  had  the  smaller  normalized  rms  error  (Figure 
8) .  The  difference  in  the  normalized  rms  error  between  T  = 
12-hr  3rd  pass  and  T  =  15-hr  3rd  pass,  was  almost  negligible. 

The  lower  model  levels  on  all  four  integration  times, 
regardless  of  the  number  of  DFI  passes,  show  an  increase  in 
the  normalized  rms  error  relative  to  the  interior  of  the 
model.  This  increase  of  error  is  not  a  result  of  a  difference 
in  the  magnitude  of  the  current,  but  rather  a  difference  in 
phase.  This  phase  error  is  illustrated  by  comparing  the  zero 
lines  of  the  vertical  velocity  profile  from  the  control  run  at 
day  10  ,to  the  corresponding  zero  lines  in  the  12-hr  3rd  pass 
profile  (Figure  9).  A  reason  for  this  error  is  that  the 
10-day  control  run  had  a  bottom  friction  term,  while  the  DFI 
procedure  did  not.  In  order  for  DFI  to  remain  adiabatic 
(friction  cannot  be  used  in  a  backward  integration)  the 
frictional  component  could  not  be  allowed  in  its  makeup. 
Therefore,  the  zero  line  for  the  12-hr  3rd  pass  DFI  solution 
is  ahead  of  that  in  the  10-day  control  run.  This  hypothesis, 
concerning  the  friction  term,  could  be  tested  by  rerunning  the 
10-day  forward  integration  without  the  friction  term,  and  then 
compare  the  vertical  velocities  (DFI  and  control)  at  the  lower 
model  levels  to  see  if  the  phase  error  had  been  eliminated. 
However,  this  is  beyond  the  scope  of  this  thesis,  and  is  left 
to  future  studies. 
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Figure  7.  Normalized  RMS  Error,  (a)  T  =  12-hr.  (b)  T  -  15-hr. 
Solid  line,  Ipaas;  dot  dash  line,  2pa3s;  dot  line,  3pass. 
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Figure  9.  Comparison  of  Control  and  DFI  w  Vertical  Profiles, 
(a)  Control  (m  d'^)  .  (b)  DFI,  T  =  12-hr  Spaas  (m  d~^)  . 


The  rms  of  the  difference  between  the  vertical 
velocities,  i.e.  without  normalization,  shows  the  error  for 
the  different  integration  times  and  the  number  of  DFI  passes 
(Figures  10,  and  11) .  The  error  decreased  with  the  number  of 
passes  and  with  the  longer  integration  times.  The  error  is 
large  where  the  vertical  velocity  is  large  and  it  has  a  sharp 
decrease  at  level  5  (250  m)  where  there  is  a  rapid  phase 
change  in  the  vertical  corresponding  to  the  bottom  of  the 
vertical  exponential  structure  (Figure  12) .  This  is 
contrasted  with  the  normalized  rms  error  that  has  a  minimum 
value  at  the  large  vertical  velocities,  but  a  sharp  increase 
at  the  vertical  phase  change  at  level  5  (Figure  13). 

From  the  results  of  the  different  integration  times  and 
the  different  number  of  repeated  DFI  passes,  it  is  clear  that 
the  12-hr  and  15-hr  integration  times  reproduce  the  control 
the  best.  However,  considering  the  increased  computational 
time  required  for  the  15-hr  DFI  integration,  it  does  not 
return  a  normalized  rms  error  that  is  significantly  better 
than  the  12-hr  DFI  integration.  The  results  also  show  that 
the  rms  of  the  difference  between  vertical  velocities  is  small 
for  the  12-hr  3rd  pass  DFI  method,  at  all  model  levels  (Figure 
14. a),  and  the  normalized  rms  error  is  under  9%  for  all 
levels,  except  the  lowest  six  (2250  -  4000  m) (Figure  14. b,  and 
Table  I) .  This  increased  error  at  the  lower  six  levels  is 
hypothesized  to  be  due  to  the  removal  of  the  friction  term  in 
the  DFI  integration.  The  end  result  is  that  the  errors  are 
acceptably  small,  and  that  the  12-hr  3rd  pass  DFI  method 
returns  the  original  values  with  highly  useful  accuracy  as 
shown  in  Figures  9,  and  15. 

As  stated  earlier,  the  DFI  procedure  was  tested  using  two 
different  initial  velocities.  The  first  series  of  tests  set 
the  horizontal  velocity  to  the  geostrophic  velocities  confuted 
from  the  density  field  at  the  10th  day  of  the  run.  This  has 
been  reviewed  above.  The  second  series  of  tests  set  the 
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Figure  10.  RMS  of  the  Difference  (m  d'^)  .  (al  T  -  6-hr.  (bj 
T  “  9-hr.  Solid  line,  Ipass;  dot  dash  line,  2pass;  dot  line, 
Spass. 
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Figure  12.  (a)  RMS  of  the  Difference  combined  plot.  Solid 
line,  6-hr  Ipass;  dot  dash  line,  9-hr  2pass;  dot  line,  12-hr 
3passf  dash  line,  15-hr  3pass.  (b)  W  vertical  profile  (m  d"‘) . 
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Figure  13.  (a)  Normalized  RMS  Error  combined  plot.  Solid  line, 
6-hr  Ipass;  dot  dash  line,  9-hr  2pa33;  dot  line,  12-hr  Spass; 
dash  line,  15-hr  Spass.  (b)  W  vertical  profile  (m  d'^)  . 


Figure  14.  DFI  Verification  Using  T  =  12-hr  3passes.  (a)  RMS 
of  the  difference  (m  d“^)  .  (b)  Normalized  RMS  Error. 
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Rossby  Wave: 


T  =  12-hr  3pass 


Model 

Level  Depth  (m) 

RMS  DIFFERENCE 

RMS(w  -  w,) 

NRE 

RMS(w  -  w,)/RMS{w) 

1 

10 

5.7684199  E-04 

5.3825520  E-02 

2 

40 

1.5399987  E-03 

4.6124276  E-02 

3 

90 

2.0126472  E-03 

3.9678775  E-02 

4 

160 

1.2841029  E-03 

2.4625398  E-02 

5 

250 

9.1466750  E-04 

2.4539951  E-02 

6 

;  360 

8.6603622  E-04 

3.0977732  E-02 

7 

490 

2.4951532  E-04 

6.3547445  E-03 

8 

640 

2.1893358  E-04 

4.8328936  E-03 

9 

1  810 

5.7628244  E-04 

1.3992130  E-02 

10  1000 

2.7732618  E-04 

8.6974138  E-03 

11  1210 

7.0427696  E-04 

3.2166429  E-02 

12  1440 

3.9193392  E-04 

2.8955661  E-02 

13  1690 

6.1763701  E-04 

8.2097113  E-02 

14  1960 

6.9397717  E-04 

0.1856848 

15  2250 

2.5081812  E-04 

0.1379997 

16  2560 

6.7364780  E-04 

0.5211930 

17  2890 

1.0700605  E-03 

0.8995180 

18  3240 

1.0504758  E-03 

1.102790 

19  3610 

6.4756384  E-04 

1.198437 

Table  I.  Error  Between  Control  and  DFI  w  Velocities. 
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Figure  15.  Comparison  of  Control  and  DFI  w  Horizontal 
Profiles,  90  m.  (a)  Control  (m  d’^) .  (b)  DFI,  T  -  12-hr  Spass 
(m  d’M  . 


initial  horizontal  velocities  to  zero,  before  running  the  DFI 
on  the  density  field.  This  was  done  to  examine  the  case  of 
using  a  very  poor  "guess"  for  the  initial  horizontal 
velocities,  and  for  theoretical  curiosity.  From  using  the 
geostrophic  velocity  as  the  initial  velocity,  it  has  been 
demonstrated  above  that  the  12-'hr  time  span  and  3rd  pass 
method  returned  the  most  accurate  velocities.  This  same  time 
span,  and  number  of  repeated  passes,  was  used  in  the  second 
DFI  test  in  which  the  initial  currents  were  set  to  zero. 

The  following  equation  shows  the  relationship  between  the 
final  filtered  (DFI)  velocity  and  the  geostrophic  velocity 
determined  by  the  initial  density  field; 

V,  -  Viy  -  (  V,  -  V,,)  +  ( V,,  -  Vjj)  (4) 

Here  V,  is  the  filtered  (dynamically  balanced)  velocity  from 
DFI,  and  is  the  geostrophic  velocity  computed  from  the 

initial  density  field  (day  10  of  the  10-day  integration) .  V,, 

is  the  (filtered)  geostrophic  velocity  computed  from  the 
initialized  density  field  using  the  DFI  method.  The  left  side 
of  Equation  (4)  gives  the  total  difference  (dynamical 
adjustment)  between  the  final  (balanced)  currents  and 
geostrophic  currents  determined  by  the  initial  mass  field. 
The  first  pair  of  terms  on  the  right  hand  aide  measure  the 
extent  to  which  the  initialized  currents  are  ageostrophic, 
while  the  last  pair  of  terms  measure  the  extent  to  which  the 
density  has  changed  during  the  initialization  process.  Which 
of  these  two  pair  of  terms  are  larger  will  depend,  obviously, 
on  the  initial  currents  used  in  the  DFI  procedure. 

The  results  (not  shown)  indicate  that  the  difference 
between  the  initial  and  final  density  fields  (right  hand  side, 
last  two  te3ntis)  was  larger  than  the  difference  due  to  the 
geostrophic  approximation  (right  hand  side,  first  two  terms) . 
This  indicates  that  the  initialized  (DFI)  density  field  was 
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different  from  the  density  field  of  the  control  run,  and  the 
initialized  currents  were  in  near  geostrophic  balance  with  the 
initialized  density  field.  That  is,  the  density  was  changed 
by  the  initialization  procedure.  This  could  be  explained  in 
terms  of  the  scales  being  examined  (length  scale  not  much 
larger  than  the  Rossby  radius),  and  that  the  disturbance  in 
the  density  field  had  flattened  out  during  initialization. 
This  flattening  could  be  the  result  of  no  support  of  the 
density  disturbance  by  the  initial  currents,  which  were  zero. 
This  particular  result  is  scale  dependent,  and  it  has  little 
practical  application  since  geostrophic  currents  can  (and 
should)  always  be  used  for  the  initial  currents  on  these 
scales.  This  issue  will  not  be  investigated  further  in  this 
study. 

B.  BABOCLDIIC  UNSTABLE  MEAN  FLOW  -  RUN  TO  FINITE  AMPLITUDE 

1 .  Setiqp 

The  baroclinic  unstable  mean  flow  was  setup  in  the  same 
manner  as  the  linear  Rossby  wave  run.  The  model  was 
initialized  using  DFI  with  a  12-hr  forward  and  backward 
integration  that  was  repeated  3  times.  The  initial  density 
field  consisted  of  a  disturbance  superimposed  on  a  baroclinic 
mean  zonal  jet.  The  disturbance  was  given  by  the  following 
equation: 


diet  - 


i  IF  ^ 

T' 50*— ♦sin(2*ii*— *2*n*znd) 
n  xm  2 


(5) 


In  Equation  (5),  the  summation  extends  from  n  =*  1  to  N  =  14, 
which  was  the  maximum  zonal  wavenumber  used.  Here  xm  »  112*ax 
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and  ym  =  56*Ay/  which  are  the  x-  and  y-  scales  of  the  domain. 
This  equation  represents  a  density  disturbance  that  is  summed 
over  14  different  wavelengths  in  the  zonal  direction.  The 
"rnd”  is  a  random  number  which  will  give  a  different  phase  to 
each  wave.  The  longest  wavelength  in  the  x~direction  is  xm  (n 
-  1),  and  the  shortest  is  8ax  (n  ■=  14)  .  The  amplitude  of  each 
wave  is  inversely  proportional  to  its  wavenumber,  so  that  the 
geostrophic  velocity  associated  with  that  wavenumber,  which  is 
proportional  to  the  x-derivative  of  the  density  field,  will 
have  an  equal  amplitude  for  each  wavelength.  The  end  result 
is  a  sinusoidal  form  of  baroclinic  white  noise  with  a  random 
phase  adjustment,  that  contains  wavenumbers  from  1  to  14. 
This  allowed  the  instability  to  select  the  most  unstable 
wavenumber.  The  form  of  the  density  disturbance  in  the  y- 
direction  is  the  hyperbolic  secant  squared  of  {(y-yd)/yl) 
which  has  a  maximum  in  the  center  of  the  domain  and  has  a 
length  scale  in  the  y-direction  of  yl.  This  length  scale  yl, 
was  set  to  25  Ion  with  the  result  that  the  magnitude  of  the 
disturbance  drops  to  10%  of  the  maximum  in  just  10  grid  points 
(40  Ion).  The  disturbance  in  Equation  (5)  was  multiplied  by  50 
to  ensure  enough  initial  energy  was  provided  for  the 
instabilities.  One  of  the  differences  between  this  simulation 
and  the  Rossby  wave  simulation,  is  the  shallower  vertical 
profile  used  here.  The  vertical  profile  of  the  density 
disturbance  was  constant  from  the  surface  to  a  depth  of  12.5 
m,  at  25  m  it  began  a  smooth  transition  into  an  exponential 
profile.  The  exponential  profile  had  a  scale  depth  of  100  m, 
wl^ich  was  140  m  shallower  than  the  Rossby  wave  run.  Another 
difference  alluded  to  earlier  was  the  model's  domain  was 
doubled  in  the  x-direction  for  this  test  (from  56  to  112  grid 
points) .  This  increase  of  length  was  to  allow  more  room  in 
which  the  different  wavelengths  could  develop.  With  ax  =  Ay 
=  4  Ion,  the  east  to  west  distance  became  448  )cm,  while  the 
north  to  south  distance  remain  at  224  km. 
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The  mean  zonal  shear  flow  that  the  density  disturbance 
was  superimposed,  was  produced  from  the  same  density  structure 
(hyperbolic  tangent)  as  in  the  Rossby  wave  run.  However,  the 
width  of  the  current  (yl)  was  reduced  from  the  1000  )an  used 
with  the  Rossby  wave,  to  25  km.  This  tightened  the  baroclinic 
flow  from  which  the  instabilities  would  grow.  The  vertical 
profile  of  the  mean  flow  (i.e.  the  axis  of  the  jet)  was  also 
shifted  toward  the  south  with  depth,  which  resulted  in  the 
mean  flow  being  centered  slightly  south  of  the  center  of  the 
y“axis.  The  southward  tilt  in  the  zonal  flow  increased  the 
instability  in  the  model  by  increasing  the  vertical  shear. 
The  increased  shear  increased  the  instabilities  in  the  model 
and  helped  the  disturbances  grow  quickly.  Besides  the  tilt, 
the  zonal  shear  flow  was  also  multiplied  by  the  same  shallow 
z-dependent  function,  as  the  density  disturbance  (constant  in 
z  above  12.5  m,  at  25  m  a  transition  to  an  exponential 
function  with  a  100  m  scale  depth) . 

The  model  was  integrated  forward  40  days  after  the  DFI 
initialization  in  order  to  produce  the  control  run.  This  was 
to  allow  the  instabilities  to  grow  and  develop  to  a  sufficient 
magnitude.  The  instantaneous  density  field  was  extracted  from 
the  40th  day,  and  used  in  the  following  DFI  tests  in  an 
attempt  to  recover  the  control  currents.  The  geostrophic 
velocities  calculated  from  the  day  40  density  structure  were 
used  as  initial  horizontal  velocities  for  the  DFI.  The  DFI 
tests  were  performed  with  three  different  integration  times  (T 
=  9,  12,  and  15  hours) ,  and  up  to  three  repeated 
initialization  passes,  for  each  of  the  integration  times. 

2 .  Control  Slaulatlon 

The  baroclinic  unstable  flow  simulation  resulted  in  a 
sinusoidal  distribution  (zonal  wavenumber  2)  of  divergence, 
vorticity  and  a  tilted  zonal  distribution  for  the  w  velocity 
at  day  40  as  shown  in  Figures  16,  and  17. a.  The  ageostrophic 
vector  field  (Figure  17. b)  shows  a  series  of  cells  (rotational 
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Figure  1«.  Unstable  Baroclinlc  Flow,  (a)  Divergence (contour 
0.8)  (lO-’s-^) .  (b)  Vorticity{ contour  20.0)  {10-''3-M  •  Negative 

values  are  dashed.  Depth  65  m. 


Figure  17.  Unstable  Baroclinic  Flow,  (a)  W  horizontal  profile 
(m  d”*),  90  m.  (b)  Ageostrophic  vector  field,  €5  m. 
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MAXIMUM  VECTOR 


flow)  which  creates  the  areas  of  maximum  divergence  that 
produces  the  corresponding  areas  of  maximum  upward  vertical 
velocity.  The  rotational  flow  also  causes  the  areas  of 
maximum  convergence  and  their  corresponding  areas  of  maximum 
downward  vertical  velocity. 

Density,  w  velocity,  and  the  v  velocity  field  all  show 
the  sinusoidal  pattern  with  the  dominant  zonal  wavenumber  2 
(Figures  18,  19. a,  and  19. b) .  The  u  velocity  pattern  has  a 
split  sinusoidal  form  with  the  isotachs  concentrating  in  or 
just  ahead  of  the  troughs  and  expanding  in  or  just  ahead  of 
the  ridges  in  the  parallel  flows  (Figure  19. c) .  The  resulting 
w  field  is  in  quasi-geostrophic  balance  (v  and  w  in  phase), 
with  rising  motion  ahead  of  the  troughs  and  sinking  motion 
ahead  of  the  ridges  as  seen  in  Figures  19. a,  and  19.b. 


3.  DFI  Verification 

Similar  to  the  Rossby  wave  tests,  the  DFI  method  was  used 
here  to  try  to  recover  the  currents  at  day  40  of  the  control 
run,  using  only  the  (control)  density  field  from  that  time. 
Utilizing  the  normalized  rms  error  (Equation  (3))  the  accuracy 
of  the  DFI  solution  was  examined.  The  DFI  was  examined  with 
regard  to  the  different  time  spans,  and  to  the  different 
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number  of  passes  of  integration. 

The  first  verification  was  made  for  just  one  DFI  pass  for 
each  of  the  three  integration  times  (T  =9,  12,  and  15  hrs) . 
For  all  time  spans,  the  control’s  u  and  v  velocity  fields 
were  almost  completely  recovered  by  the  first  DFI  pass. 
However,  the  w  field  required  the  integration  time  to  be 
lengthened,  from  T=9toT*15  hrs,  in  order  for  more  of  it 
to  be  returned  (Figure  20) .  As  shown  in  Figure  20,  there  is 
a  dramatic  decrease  in  the  normalized  rms  error  when  T  is 
increased  from  9-hr  to  12-hr.  There  is  a  smaller  increase  in 
accuracy  between  the  12-hr  and  15-hr  integration.  As  in  the 
Rossby  wave  case,  the  error  remains  significant  (NRE  -  0.4) 
near  the  bottom.  However,  in  this  case  it  is  also  higher  near 
the  surface.  The  interpretation  of  this  error  is  given  below. 

The  DFI  procedure  was  then  repeated  up  to  three  times, 
and  the  results  were  compared  for  each  integration  time. 
These  results  show  that  for  the  three  different  integration 
times  (T  =  9,  12,  and  15  hrs),  the  normalized  rms  error 
between  the  w's,  decreased  with  the  number  of  DFI  passes 
(Figures  21  and  22)  .  The  least  improvement  with  regard  to  the 
normalized  rms  error  was  found  in  the  second  and  third  passes. 
Notable  are  the  second  and  third  passes  of  the  15-hr 
integration  (Figure  22) .  This  indicate  that  after  the  second 
DFI  pass  almost  all  of  the  control  vertical  velocity  has  been 
recovered.  Figure  23  shows  that  for  the  different  time  spans 
at  the  3rd  pass,  the  longer  integration  times  had  the  smaller 
normalized  rms  error.  The  difference  between  T  =  12-hr  3rd 
pass  and  T  -  15-hr  3rd  pass  with  regard  to  the  normalized  rms 
error  was  insignificant. 

The  lower  model  levels,  for  all  four  integration  times, 
show  an  increase  in  the  magnitude  of  the  normalized  rms  error. 
This  increase  of  error  is  due  to  a  combination  of  effects. 
First,  there  is  a  slight  phase  error  in  the  vertical 
velocities.  By  comparing  the  zero  lines  of  the  control  and 
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Figure  20.  Normalized  RMS  Error.  First  pass  with  different 
integration  times.  Solid  line,  9-hr;  dot  dash  line,  12-hr;  dot 
line,  15-hr. 
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Ufels  ItTels 


Figure  22.  Normalized  RMS  Error,  (a)  T  -=  15-hr.  Solid  line, 
Ipass;  dot  dash  line,  2pas3;  dot  line,  3pass. 
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Flguro  23.  Normalized  RMS  Error.  Third  pass  with  different 
integration  times.  Solid  line,  9-hr;  dot  dash  line,  12-hr;  dot 
line,  15-hr. 

DFI  vertical  velocity,  in  Figure  24,  it  is  seen  that  near  the 
western  edge  of  the  domain  the  DFI  zero  line  is  slightly  ahead 
of  the  control,  but  in  the  middle  of  the  domain  they  are  the 
same,  ■  then  the  control  moves  ahead  of  the  DFI  at  the  eastern 
edge.  One  would  have  thought  because  the  40-day  control  run 
had  a  bottom  friction  term,  while  the  DFI  procedure  did  not, 
that  the  DFI  zero  line  would  consistently  be  ahead  of  the 
control's  zero  line  as  in  the  Rossby  wave  case.  Instead,  only 
part  of  the  normalized  rms  error  is  due  to  a  phase  error. 

The  second  part  of  this  normalized  rms  error  can  be  found 
by  looking  again  at  Figure  24.  There  it  is  seen  that  the  DFI 
run  has  positive  contours  going  slightly  deeper  than  the 
contours  for  the  control  at  the  lower  levels.  The  negative 
contours  also  extend  deeper.  This  would  indicate  a  slight 
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magnitude  error  in  the  recovery  of  the  currents  at  these 
depths  with  the  DFI  solution  somewhat  stronger.  This  is  also 
seen  by  the  larger  central  values  (amplitudes)  for  the  DFI 
solution.  A  rapid  increase  in  the  percent  of  error  is  also 
observed  near  the  surface.  It  appears  here  too,  that  the 
positive  contours  for  the  DFI  are  just  a  bit  larger  than  the 
control's.  These  larger  errors,  that  occur  in  both  the 
shallower  and  deeper  parts  of  the  domain,  are  located  where 
the  model's  resolution  of  the  phenomenon  is  poorest.  For  the 
error  near  the  surface,  there  is  a  shallow  vertical  scale  that 
only  extends  over  1  or  2  model  levels.  In  the  deeper  ocean, 
the  vertical  levels  are  grouped  farther  and  farther  apart. 
This  hypothesis  could  be  addressed  by  increasing  the  vertical 
resolution.  The  simulation  could  be  done  again  with  the 
results  compared  between  the  control  and  DFI  to  see  if  the 
normalized  rms  error  had  decreased.  This  was  not  done  and 
left  for  future  study.  Another  possible  explanation  is  the 
fact  that  the  DFI  solution  is  adiabatic  (no  friction),  so  the 
amplitude  of  the  disturbances  are  not  damped  as  in  the  control 
run.  Comparing  the  central  values  of  Figures  24  show  the 
larger  magnitudes  of  the  DFI  currents.  The  actual  test  of 
this  premise  would  require  a  frictionless  control  run,  which 
was  also  left  for  a  future  study. 

Though  the  normalized  rms  error  increased  at  the  surface 
and  at  the  lower  7  levels,  the  absolute  error  at  these 
locations  is  small  for  all  time  spans.  The  rms  of  the 
difference  between  the  vertical  velocities  (i.e.  not 
normalized) ,  shows  the  absolute  error  for  the  different 
integration  times  and  the  different  number  of  DFI  passes 
(Figures  25,  and  26) .  The  absolute  error  decreased  with  the 
number  of  passes  and  with  the  longer  integration  times.  The 
absolute  error  is  large  where  the  vertical  velocity  is  large 
at  level  3  (90  m)  (Figure  27).  This  behavior  is  contrasted 
with  that  of  the  normalize  rms  error  that  has  a  minim\am  value 
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Figure  25.  RMS  of  the  Difference  (m  d"^)  .  (a)  T  -  9-hr.  (b)  T 
12-hr.  Solid  line,  Ipass;  dot  dash  line,  2pass;  dot  line, 
3pass. 
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Flyure  27.  (a)  RMS  of  the  Difference  combined  plot  (m  d"^)  . 
Solid  line,  9--hr  Ipass;  dot  dash  line,  12-hr  2pa3s;  dot  line, 
15-hr  3pas3.  (b)  W  vertical  profile  (m  d“‘)  . 
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where  the  vertical  velocities  are  large,  but  increases  where 
the  values  of  w  are  small  (Figure  28). 

Li)ce  the  findings  from  the  Rossby  wave  tests,  the  results 
show  that  the  12-hr  and  15-hr  integration  times  reproduce  the 
control's  currents  the  best.  However,  the  increased 
computational  time  needed  for  the  15-hr  DFI  integration,  does 
not  produce  a  significantly  smaller  percent  of  error. 
Therefore,  the  DFI  procedure  with  the  12-hr  integration  time 
and  3  passes  is  considered  to  give  the  best  practical  results. 

The  results  show  that  the  rms  of  the  difference  between 
vertical  velocities  is  small  for  the  12— hr  3rd  pass  DFI 
method,  at  all  model  levels  (Figure  29. a),  and  the  normalized 
rms  error  is  under  15%  for  all  levels,  except  at  the  surface 
and  the  lowest  five  levels  (2560  -  4000  m) (Figure  29. b,  and 
Table  II).  This  increased  error  near  the  top  and  bottom  is 
most  likely  due  to  the  lack  of  resolution  in  these  regions  and 
the  adiabatic  assumption  of  the  DFI  process.  The  conclusion 
is  that  the  errors  are  acceptably  small,  and  that  the  12-hr 
3rd  pass  DFI  method  returns  the  original  values  with  highly 
useful  accuracy  as  shown  in  Figures  24,  and  30. 
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Figure  28.  (a)  Normalized  RMS  Error  combined  plot.  Solid  line, 
9-hr  Ipass;  dot  dash  line,  12-hr  2pass;  dot  line,  15-hr  3pass. 
(b)  W  vertical  profile  (m  d"*). 
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Imls  l««l» 


Baroclinic  Unstable  Flow: 


12-hr  3pass 


Model 


RMS  DIFFERENCE 


Level  Depth (m)  RMS(w  -  w,) 


0.1553140 


0.2500602 


0.2185996 


0.1549404 


0.1165265 


0.1109206 


0.1107548 


0.1167470 


.01180682 


0.1211561 


0.1237135 


0.1264377 


0.1307218 


0.1338995 


0.1335814 


0.1308055 


0.1236864 


0.1082318 


7.0083916  E-02 


RMS(w  -  w,)  /  RMS(w) 


1000 


1210 


1440 


1690 


1960 


2250 


2560 


2890 


3240 


3610 


0.1408043 


9.2744939  E-02 


6.7110352  E-02 


5.7649828  E-02 


6.1324663  E-02 


6.5937974  E-02 


7.3884249  E-02 


7.9604834  E-02 


8.7748334  E-02 


9.7322747  E-02 


0.1095053 


0.1267268 


0.1484055 


0.1741401 


0.2089019 


0.2580079 


0.3314812 


0.4196966 


Table  II.  Error  Between  Control  and  DFI  w  Velocities. 
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Figure  30.  Comparison  of  Control  and  DFI  w  Horizontal 
Profiles,  90  m.  (a)  Control  (m  d'")  .  (b)  DFI,  T  -  12-hr  Spasses 
{m  d"^)  . 
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III.  DISCUSSION  AMD  CONCLUSIONS 


The  DFI  method  was  tested  to  see  how  accurately  it  could 
recover  the  three  dimensional  circulation  from  the  initial 
density  field  and  its  geostrophic  velocities  alone.  Two  types 
of  disturbances  were  used,  a  linear  Rossby  wave,  and  a  finite 
amplitude  disturbance  growing  on  an  unstable  baroclinic  flow. 
Four  different  integration  times  were  used  in  the  DFI 
procedure  with  the  Rossby  wave  and  three  with  the  baroclinic 
unstable  flow,  all  of  which  where  repeated  up  to  three  times. 
The  absolute  error  and  the  percent  error  for  each  time  span 
and  pass  were  calculated.  These  "measures  of  skill"  were  then 
compared  to  determine  which  best  diagnosed  the  currents  from 
the  control  run.  Attention  was  focused  on  the  vertical 
velocity  because  of  its  large  signal  and  its  dynamical 
importance . 

In  both  test  cases  the  12-hr  integration  time  with  3 
passes  produced  the  best  practical  results.  Though  there 
continued  to  be  increases  of  accuracy  at  longer  integration 
times,  these  were  small  enough  to  be  confident  that  the  12-hr 
integration  had  recover  the  majority  of  the  current.  There 
were  certain  preferred  areas  that  exhibited  larger  normalized 
rms  errors  (NEIE  ~  0.3)  in  both  cases. 

The  Rossby  wave  had  increased  normalized  rms  errors  at 
the  lower  levels.  This  was  thought  to  be  the  result  of  the 
friction  term  used  in  the  control  run,  but  not  in  the  DFI. 
The  DFI  run  was  constrained  to  remain  adiabatic  because  it 
could  not  include  such  a  term  in  its  backward  integration. 
This  produced  a  phase  error  between  the  control  and  the  DFI, 
where  the  DFI  was  60-70  degrees  ahead  of  the  control.  New 
work  by  Lynch  and  Huang  are  removing  the  requirement  for  the 
backward  integration  by  the  use  of  new  filters  that  would  use 
only  forward  integration  and  still  be  able  to  pull  the 
filtered  fields  back  to  t  -  o.  This  would  remove  the 
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requirement  for  an  adiabatic  process  and  allow  friction  to  be 
included  in  the  DFI  method.  The  inclusion  of  the  friction 
term  would  mean  a  more  accurate  diagnoses  of  the  actual  ocean 
environment,  where  friction  is  found  not  only  at  the 
fluid/solid  interface  (coastline  and  ocean  floor),  but  also  in 
the  upper  levels  where  the  fluid/air  interface  occur  with 
resulting  mixing.  For  the  Rossby  wave  the  DFI  returns  the 
control  currents  with  less  than  10%  error,  except  where  noted. 
If  the  friction  term  is  included,  a  10%  or  less  error  should 
be  extended  to  the  total  domain. 

Another  possible  explanation  for  the  error  that  is  found 
at  the  lower  levels  in  the  domain,  is  that  the  stability  there 
is  almost  zero  which  may  haitqper  the  geostrophic  adjustment. 
If  the  lower  levels  are  unable  to  adjust  then  no  matter  what 
fix  is  applied  to  the  model,  there  will  always  remain  some 
error.  The  baroclinic  unstable  flow  had  significant 
normalized  rms  errors  both  near  the  surface  and  at  the  lower 
levels.  These  errors  could  be  the  result  of  poor  vertical 
resolution  and/or  be  some  of  the  same  problems  due  to  the 
adiabatic  (frictionless)  assumption  in  the  DFI.  The 
relatively  poor  resolution  near  the  surface  was  induced  by  the 
shallow  vertical  structure  of  the  density  disturbance  where 
only  a  few  levels  were  able  to  resolve  it.  The  lower  level 
resolution  error  maybe  the  result  of  the  larger  spacing 
between  the  levels  with  depth.  The  error  could  also  be  the 
result  of  the  frictionless  DFI  procedures,  or  the  failure  of 
the  geostrophic  adjustment  process  in  the  deep  ocean  where  the 
stability  is  so  weak.  It  is  hoped  that  at  least  two  of  these 
errors  can  be  reduced  by  Increasing  the  vertical  resolution, 
and  by  running  the  control  without  friction.  This  could  be 
done  simply  by  reducing  the  total  model  depth  and  removing  the 
frictional  component  from  the  control  run. 

Despite  these  errors,  the  DFI  method  in  fact  returned  the 
control’s  currents,  including  the  secondary  circulation  and 
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vertical  velocity,  with  a  high  degree  of  accuracy.  What  is 
significant  is  that  the  vertical  velocity  was  accurately  found 
by  using  only  the  density  field,  and  that  it  was  done  quickly 
and  simply.  This  has  great  implications  to  forecast  models, 
because  now  all  velocities  can  be  accurately  found  for  a  given 
location  which  can  then  be  used  to  initialize  the  model.  Also 
the  recovered  velocities  have  had  almost  all  high  frequency 
noise  removed,  which  reduces  the  problems  that  gravity  waves 
have  caused  in  forecast  models  in  the  past.  By  applying  the 
DFI  to  these  models  it  will  increase  the  confidence  of  their 
results  and  increase  the  understanding  and  knowledge  of  the 
mesoscale  features  in  the  littoral  zone. 
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